# Replication File for
# "Coalition Prospects and Policy Change: An Application to the Environment"
# Authors: Mark A. Kayser and Jochen Rehmert


# load required packages
library(ggplot2);library(grid);library(stargazer);library(plm);library(lmtest);library(sandwich)

# set directory
setwd('C:/Users/rehmerjo/Dropbox/CIP_LSQ_Conference/dataverse')

# read data
dat <- read.csv("EPS_data.csv", sep = ";", dec = ",")

# NB: The file has been saved on a German machine.
# This means that any .csv file will use the ";" as 
# a delimiter and the "," as the decimal separator

# coding of control variables 
dat$kyoto <- 0
dat$kyoto[dat$country == "austria" & dat$year >= 1998] <- 1
dat$kyoto[dat$country == "czechia" & dat$year >= 1998] <- 1
dat$kyoto[dat$country == "denmark" & dat$year >= 1998] <- 1
dat$kyoto[dat$country == "estonia" & dat$year >= 1998] <- 1
dat$kyoto[dat$country == "finland" & dat$year >= 1998] <- 1
dat$kyoto[dat$country == "germany" & dat$year >= 1998] <- 1
dat$kyoto[dat$country == "greece" & dat$year >= 1998] <- 1
dat$kyoto[dat$country == "ireland" & dat$year >= 1998] <- 1
dat$kyoto[dat$country == "italy" & dat$year >= 1998] <- 1
dat$kyoto[dat$country == "japan" & dat$year >= 1998] <- 1
dat$kyoto[dat$country == "netherlands" & dat$year >= 1998] <- 1
dat$kyoto[dat$country == "new zeland" & dat$year >= 1998] <- 1
dat$kyoto[dat$country == "norway" & dat$year >= 1998] <- 1
dat$kyoto[dat$country == "poland" & dat$year >= 1998] <- 1
dat$kyoto[dat$country == "portugal" & dat$year >= 1998] <- 1
dat$kyoto[dat$country == "slovakia" & dat$year >= 1999] <- 1
dat$kyoto[dat$country == "slovenia" & dat$year >= 1998] <- 1
dat$kyoto[dat$country == "spain" & dat$year >= 1998] <- 1
dat$kyoto[dat$country == "sweden" & dat$year >= 1998] <- 1

dat$fukushima <- 0
dat$fukushima[dat$year >= 2011] <- 1

##################################################################
# Figure A1                                                      #
##################################################################

### CIP & Polls ###
# fixed-effects, clustered SE
summary(mod.1.a <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       eco_cip_gross + I(eco_polls/100) , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0 ,]))
mod.1.a.se <- coeftest(mod.1.a, vcov=vcovHC(mod.1.a,type="HC0",cluster="group"))[,2]


# lagged dependent variable + trend
summary(mod.1.b <-lm(eps ~  l_eps + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg +  eco_gov + environ_eco +
                        eco_cip_gross + I(eco_polls/100) + year  , 
                      data = dat[!is.na(dat$eco_polls) & dat$cee == 0,]))
mod.1.b.se <- coeftest(mod.1.b, vcov=vcovHC(mod.1.b,type="HC0",cluster="group"))[,2]


# first-difference
summary(mod.1.c <-plm(eps ~  cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                        eco_cip_gross +  I(eco_polls/100)  , 
                      data = dat[!is.na(dat$eco_polls) & dat$cee == 0,], model = "fd", index = c("country_id","year")))
mod.1.c.se <- coeftest(mod.1.c, vcov=vcovHC(mod.1.c,type="HC0",cluster="group"))[,2]



stargazer(list(mod.1.a, mod.1.b, mod.1.c), se = list(mod.1.a.se, mod.1.b.se, mod.1.c.se), 
          omit = c("country_id","year"),
          covariate.labels = c("EPS$_{t-1}$", "Minority Cabinet", "Kyoto Protocol","GDP Growth","Cabinet's Mean Environmental Protection",
                               "Total Greenhouse Gas Emissions/Capita", "Green Party in Government",  "Green Party's Environmental Protection", "Green Party's gross CIP",
                               "Green Party's Polls", "Intercept"),
          add.lines = list(c("Number of Countries", length(unique(mod.1.a[["model"]]$`as.factor(country_id)`)), length(unique(mod.1.a[["model"]]$`as.factor(country_id)`)), length(unique(mod.1.a[["model"]]$`as.factor(country_id)`))),
                           c("Country Fixed-Effects", "Yes", "No", "No"),
                           c("Trend", "No","Yes","No")))

##################################################################
# Figure A2                                                      #
##################################################################

# fixed-effects, clustered SE
summary(mod.1.a <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       eco_cip_gross  , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0 ,]))

# lagged dependent variable + trend
summary(mod.1.b <-lm(eps ~  l_eps + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg +  eco_gov + environ_eco +
                        eco_cip_gross + year  , 
                      data = dat[!is.na(dat$eco_polls) & dat$cee == 0,]))

# first-difference
summary(mod.1.c <-plm(eps ~  cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                        eco_cip_gross  , 
                      data = dat[!is.na(dat$eco_polls) & dat$cee == 0,], model = "fd", index = c("country_id","year")))


stargazer(list(mod.1.a, mod.1.b, mod.1.c), 
          omit = c("country_id","year"),
          covariate.labels = c("EPS$_{t-1}$", "Minority Cabinet", "Kyoto Protocol","GDP Growth","Cabinet's Mean Environmental Protection",
                               "Total Greenhouse Gas Emissions/Capita", "Green Party in Government",  "Green Party's Environmental Protection", "Green Party's gross CIP",
                                "Intercept"),
          add.lines = list(c("Number of Countries", length(unique(mod.1.a[["model"]]$`as.factor(country_id)`)), length(unique(mod.1.a[["model"]]$`as.factor(country_id)`)), length(unique(mod.1.a[["model"]]$`as.factor(country_id)`))),
                           c("Country Fixed-Effects", "Yes", "No", "No"),
                           c("Trend", "No","Yes","No")))

##################################################################
# Figure A3                                                      #
##################################################################

# fixed-effects, clustered SE
summary(mod.3.a <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       pm_cip_eco, 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0,]))

# pooled with LDV
summary(mod.3.c <-lm(eps ~ l_eps + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       pm_cip_eco  + year , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0,]))

# first-difference
summary(mod.3.b <-plm(eps ~  cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                        pm_cip_eco , 
                      data = dat[!is.na(dat$eco_polls) & dat$cee == 0,], model = "fd", index = c("country_id","year")))

# fixed-effects, clustered SE & interaction 
summary(mod.3.d <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       pm_cip_eco * eco_cip_gross, 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0,]))

# fixed-effects, clustered SE & interaction & trend
summary(mod.3.e <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       pm_cip_eco * eco_cip_gross + year, 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0,]))

stargazer(list(mod.3.a, mod.3.c, mod.3.b, mod.3.d, mod.3.e), 
          omit = c("country_id","year"),
          covariate.labels = c("EPS$_{t-1}$", "Minority Cabinet","Kyoto Protocol",
                               "GDP Growth","Cabinet's Mean Environmental Protection",
                               "Total Greenhouse Gas Emissions/Capita","Green Party in Government",
                               "Green Party's Environmental Protection",
                               "PM Party's CIP net of Green parties (PM CIP)", "Green Party's gross CIP (Green CIP)", "PM CIP $\\times$ Green CIP", "Intercept"),
          add.lines = list(c("Number of Countries","9","9","9","9","9"),
                           c("Fixed-Effects", "Yes", "No","No","Yes","Yes"),
                           c("Lagged Dependent Vairable","No","Yes","No","No","No"),
                           c("Trend", "No", "Yes", "No", "No", "Yes") ), dep.var.labels = "")

##################################################################
# Figure A4                                                      #
##################################################################

# fixed-effects, clustered SE
summary(mod.1.a <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       fukushima + eco_cip_gross  , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0 ,]))
mod.1.a.se <- coeftest(mod.1.a, vcov=vcovHC(mod.1.a,type="HC0",cluster="group"))[,2]


# lagged dependent variable + trend
summary(mod.1.b <-lm(eps ~  l_eps + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg +  eco_gov + environ_eco +
                        fukushima + eco_cip_gross + year  , 
                      data = dat[!is.na(dat$eco_polls) & dat$cee == 0,]))
mod.1.b.se <- coeftest(mod.1.b, vcov=vcovHC(mod.1.b,type="HC0",cluster="group"))[,2]


# first-difference
summary(mod.1.c <-plm(eps ~  cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                        fukushima + eco_cip_gross  , 
                      data = dat[!is.na(dat$eco_polls) & dat$cee == 0,], model = "fd", index = c("country_id","year")))
mod.1.c.se <- coeftest(mod.1.c, vcov=vcovHC(mod.1.c,type="HC0",cluster="group"))[,2]


stargazer(list(mod.1.a, mod.1.b, mod.1.c), 
          se = list(mod.1.a.se, mod.1.b.se, mod.1.c.se),
          omit = c("country_id","year"),
          covariate.labels = c("EPS$_{t-1}$", "Minority Cabinet", "Kyoto Protocol","GDP Growth","Cabinet's Mean Environmental Protection",
                               "Total Greenhouse Gas Emissions/Capita", "Green Party in Government",  "Green Party's Environmental Protection", "Post-Fukushima","Green Party's gross CIP",
                                "Intercept"),
          add.lines = list(c("Number of Countries", length(unique(mod.1.a[["model"]]$`as.factor(country_id)`)), length(unique(mod.1.a[["model"]]$`as.factor(country_id)`)), length(unique(mod.1.a[["model"]]$`as.factor(country_id)`))),
                           c("Country Fixed-Effects", "Yes", "No", "No"),
                           c("Trend", "No","Yes","No")))


##################################################################
# Figure A5                                                      #
##################################################################

# fixed-effects, clustered SE
summary(mod.3.a <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       fukushima +pm_cip_eco, 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0,]))
mod.3.a.se <- coeftest(mod.3.a, vcov=vcovHC(mod.3.a,type="HC0",cluster="group"))[,2]


# pooled with LDV
summary(mod.3.c <-lm(eps ~ l_eps + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       fukushima +pm_cip_eco  + year , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0,]))
mod.3.c.se <- coeftest(mod.3.c, vcov=vcovHC(mod.3.c,type="HC0",cluster="group"))[,2]

# first-difference
summary(mod.3.b <-plm(eps ~  cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                        fukushima +pm_cip_eco , 
                      data = dat[!is.na(dat$eco_polls) & dat$cee == 0,], model = "fd", index = c("country_id","year")))
mod.3.b.se <- coeftest(mod.3.b, vcov=vcovHC(mod.3.b,type="HC0",cluster="group"))[,2]

# fixed-effects, clustered SE & interaction 
summary(mod.3.d <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       fukushima +pm_cip_eco * eco_cip_gross, 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0,]))
mod.3.d.se <- coeftest(mod.3.d, vcov=vcovHC(mod.3.d,type="HC0",cluster="group"))[,2]

# fixed-effects, clustered SE & interaction & trend
summary(mod.3.e <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       fukushima +pm_cip_eco * eco_cip_gross + year, 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0,]))
mod.3.e.se <- coeftest(mod.3.e, vcov=vcovHC(mod.3.e,type="HC0",cluster="group"))[,2]

stargazer(list(mod.3.a, mod.3.c, mod.3.b, mod.3.d, mod.3.e), 
          se = list(mod.3.a.se, mod.3.c.se, mod.3.b.se, mod.3.d.se, mod.3.e.se),
          omit = c("country_id","year"),
          covariate.labels = c("EPS$_{t-1}$", "Minority Cabinet","Kyoto Protocol",
                               "GDP Growth","Cabinet's Mean Environmental Protection",
                               "Total Greenhouse Gas Emissions/Capita","Green Party in Government",
                               "Green Party's Environmental Protection", "Post-Fukushima",
                               "PM Party's CIP net of Green parties (PM CIP)", "Green Party's gross CIP (Green CIP)", "PM CIP $\\times$ Green CIP", "Intercept"),
          add.lines = list(c("Number of Countries","9","9","9","9","9"),
                           c("Fixed-Effects", "Yes", "No","No","Yes","Yes"),
                           c("Lagged Dependent Vairable","No","Yes","No","No","No"),
                           c("Trend", "No", "Yes", "No", "No", "Yes") ), dep.var.labels = "")

##################################################################
# Figure A6                                                      #
##################################################################


# fixed-effects, clustered SE
summary(mod.1.a <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       environ_worry_interpolated + eco_cip_gross  , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0 ,]))
mod.1.a.se <- coeftest(mod.1.a, vcov=vcovHC(mod.1.a,type="HC0",cluster="group"))[,2]


# lagged dependent variable + trend
summary(mod.1.b <-lm(eps ~  l_eps + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg +  eco_gov + environ_eco +
                        environ_worry_interpolated + eco_cip_gross + year  , 
                      data = dat[!is.na(dat$eco_polls) & dat$cee == 0,]))
mod.1.b.se <- coeftest(mod.1.b, vcov=vcovHC(mod.1.b,type="HC0",cluster="group"))[,2]


# first-difference
summary(mod.1.c <-plm(eps ~  cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                        environ_worry_interpolated + eco_cip_gross  , 
                      data = dat[!is.na(dat$eco_polls) & dat$cee == 0,], model = "fd", index = c("country_id","year")))
mod.1.c.se <- coeftest(mod.1.c, vcov=vcovHC(mod.1.c,type="HC0",cluster="group"))[,2]


stargazer(list(mod.1.a, mod.1.b, mod.1.c), 
          se = list(mod.1.a.se, mod.1.b.se, mod.1.c.se),
          omit = c("country_id","year"),
          covariate.labels = c("EPS$_{t-1}$", "Minority Cabinet", "Kyoto Protocol","GDP Growth","Cabinet's Mean Environmental Protection",
                               "Total Greenhouse Gas Emissions/Capita", "Green Party in Government",  "Green Party's Environmental Protection", "Share in Population worried about Environment","Green Party's gross CIP",
                                "Intercept"),
          add.lines = list(c("Number of Countries", length(unique(mod.1.a[["model"]]$`as.factor(country_id)`)), length(unique(mod.1.a[["model"]]$`as.factor(country_id)`)), length(unique(mod.1.a[["model"]]$`as.factor(country_id)`))),
                           c("Country Fixed-Effects", "Yes", "No", "No"),
                           c("Trend", "No","Yes","No")))

##################################################################
# Figure A7                                                      #
##################################################################


# fixed-effects, clustered SE
summary(mod.3.a <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       environ_worry_interpolated +pm_cip_eco, 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0,]))
mod.3.a.se <- coeftest(mod.3.a, vcov=vcovHC(mod.3.a,type="HC0",cluster="group"))[,2]


# pooled with LDV
summary(mod.3.c <-lm(eps ~ l_eps + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       environ_worry_interpolated +pm_cip_eco  + year , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0,]))
mod.3.c.se <- coeftest(mod.3.c, vcov=vcovHC(mod.3.c,type="HC0",cluster="group"))[,2]

# first-difference
summary(mod.3.b <-plm(eps ~  cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                        environ_worry_interpolated +pm_cip_eco , 
                      data = dat[!is.na(dat$eco_polls) & dat$cee == 0,], model = "fd", index = c("country_id","year")))
mod.3.b.se <- coeftest(mod.3.b, vcov=vcovHC(mod.3.b,type="HC0",cluster="group"))[,2]

# fixed-effects, clustered SE & interaction 
summary(mod.3.d <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       environ_worry_interpolated +pm_cip_eco * eco_cip_gross, 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0,]))
mod.3.d.se <- coeftest(mod.3.d, vcov=vcovHC(mod.3.d,type="HC0",cluster="group"))[,2]

# fixed-effects, clustered SE & interaction & trend
summary(mod.3.e <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       environ_worry_interpolated +pm_cip_eco * eco_cip_gross + year, 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0,]))
mod.3.e.se <- coeftest(mod.3.e, vcov=vcovHC(mod.3.e,type="HC0",cluster="group"))[,2]

stargazer(list(mod.3.a, mod.3.c, mod.3.b, mod.3.d, mod.3.e), 
          se = list(mod.3.a.se, mod.3.c.se, mod.3.b.se, mod.3.d.se, mod.3.e.se),
          omit = c("country_id","year"),
          covariate.labels = c("EPS$_{t-1}$", "Minority Cabinet","Kyoto Protocol",
                               "GDP Growth","Cabinet's Mean Environmental Protection",
                               "Total Greenhouse Gas Emissions/Capita","Green Party in Government",
                               "Green Party's Environmental Protection", "Share in Population worried about Environment",
                               "PM Party's CIP net of Green parties (PM CIP)", "Green Party's gross CIP (Green CIP)", "PM CIP $\\times$ Green CIP", "Intercept"),
          add.lines = list(c("Number of Countries","9","9","9","9","9"),
                           c("Fixed-Effects", "Yes", "No","No","Yes","Yes"),
                           c("Lagged Dependent Vairable","No","Yes","No","No","No"),
                           c("Trend", "No", "Yes", "No", "No", "Yes") ), dep.var.labels = "")

##################################################################
# Figure A8                                                      #
##################################################################

boots <- function(mod, n.sim = 1000, seed = 1104){
  out <- c()
  set.seed(seed)
  form <- as.character(as.formula(mod$call[2])[3])
  dv <- as.character(as.formula(mod$call[[2]])[2])
  form <- as.formula(paste0(dv, "~", form))
  dat.tmp <- mod[["model"]]
  colnames(dat.tmp)[which(colnames(dat.tmp)=="as.factor(country_id)")] <- "country_id"
  for(ii in 1:n.sim){
    data.rs <- dat.tmp[sample(1:nrow(dat), nrow(dat), replace = TRUE),]
    summary(mod.tmp <- lm(form, data = data.rs))
    ses <- summary(mod.tmp)$coefficients[,2]
    out <- rbind(out, ses)
    cat(ii, "\n")  
  }
return(out)
}

# fixed-effects, clustered SE
summary(mod.1.a <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       eco_cip_gross + year , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0 ,]))
mod.1.a.se <- coeftest(mod.1.a, vcov=vcovHC(mod.1.a,type="HC0",cluster="group"))[,2]
mod.1.a.se.bs <- apply( boots(mod = mod.1.a), 2, mean)
mod.1.a.no <- summary(mod.1.a)$coefficient[,2]

# with polling data
summary(mod.1.b <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       eco_cip_gross + I(eco_polls/100)+ year , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0 ,]))
mod.1.b.se <- coeftest(mod.1.b, vcov=vcovHC(mod.1.b,type="HC0",cluster="group"))[,2]

# interaction government status
summary(mod.1.c <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg  +  environ_eco +
                       eco_cip_gross * eco_gov + year , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0 ,]))
mod.1.c.se <- coeftest(mod.1.c, vcov=vcovHC(mod.1.c,type="HC0",cluster="group"))[,2]

# with Fukushima
summary(mod.1.d <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       fukushima + eco_cip_gross + year , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0 ,]))
mod.1.d.se <- coeftest(mod.1.d, vcov=vcovHC(mod.1.d,type="HC0",cluster="group"))[,2]

# with population worried
summary(mod.1.e <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       environ_worry_interpolated + eco_cip_gross + year , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0 ,]))
mod.1.e.se <- coeftest(mod.1.e, vcov=vcovHC(mod.1.e,type="HC0",cluster="group"))[,2]

# interaction time to election
summary(mod.1.f <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       time_election * eco_cip_gross + year  , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0 ,]))
mod.1.f.se <- coeftest(mod.1.f, vcov=vcovHC(mod.1.f,type="HC0",cluster="group"))[,2]




stargazer(list(mod.1.a, mod.1.a, mod.1.a, mod.1.b,mod.1.c, mod.1.d,mod.1.e,mod.1.f), 
          se = list(mod.1.a.se,mod.1.a.no, mod.1.a.se.bs, mod.1.b.se,mod.1.c.se,mod.1.d.se,mod.1.e.se,mod.1.f.se),
          omit = c("country_id","year"),
          covariate.labels = c("Minority Cabinet","Kyoto Protocol",
                               "GDP Growth","Cabinet's Mean Environmental Protection",
                               "Total Greenhouse Gas Emissions/Capita","Green Party in Government",
                               "Green Party's Environmental Protection", "Post-Fukushima", 
                               "Share in Population worried about Environment", "Time to next Election (in Years)",
                               "Green Party's gross CIP",  "Green Party's Polls",
                               "Green CIP $\\times$ Green Party in Government", 
                               "Green CIP $\\times$ Time to Election", 
                                "Intercept"),
          add.lines = list(c("Number of Countries","9","9","9", "9", "9", "9", "9"),
                           c("Fixed-Effects", "Yes", "Yes","Yes", "Yes", "Yes", "Yes", "Yes", "Yes"),
                           c("Trend", "Yes", "Yes","Yes", "Yes", "Yes",  "Yes","Yes", "Yes"),
                           c("SEs", "Clustered", "Not Clustered", "Bootstrapped", "Clustered", "Clustered", "Clustered", "Clustered", "Clustered")), dep.var.labels = "")


##################################################################
# Figure A9                                                      #
##################################################################

# fixed-effects, clustered SE
summary(mod.1.a <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       pm_cip_eco + year , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0 ,]))
mod.1.a.se <- coeftest(mod.1.a, vcov=vcovHC(mod.1.a,type="HC0",cluster="group"))[,2]
mod.1.a.se.bs <- apply( boots(mod = mod.1.a), 2, mean)
mod.1.a.no <- summary(mod.1.a)$coefficient[,2]

# with polling data
summary(mod.1.b <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       pm_cip_eco + eco_polls+ year , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0 ,]))
mod.1.b.se <- coeftest(mod.1.b, vcov=vcovHC(mod.1.b,type="HC0",cluster="group"))[,2]

# interaction government status
summary(mod.1.c <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       pm_cip_eco * eco_cip_gross + year , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0 ,]))
mod.1.c.se <- coeftest(mod.1.c, vcov=vcovHC(mod.1.c,type="HC0",cluster="group"))[,2]

# with Fukushima
summary(mod.1.d <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       fukushima + pm_cip_eco + year , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0 ,]))
mod.1.d.se <- coeftest(mod.1.d, vcov=vcovHC(mod.1.d,type="HC0",cluster="group"))[,2]

# with population worried
summary(mod.1.e <-lm(eps ~  as.factor(country_id) + cab_minority + kyoto+ gdp_quarter + environ_cab_mean + total_ghg + eco_gov +  environ_eco +
                       environ_worry_interpolated + pm_cip_eco + year , 
                     data = dat[!is.na(dat$eco_polls) & dat$cee == 0 ,]))
mod.1.e.se <- coeftest(mod.1.e, vcov=vcovHC(mod.1.e,type="HC0",cluster="group"))[,2]



stargazer(list(mod.1.a, mod.1.a, mod.1.a, mod.1.b,mod.1.c, mod.1.d,mod.1.e), 
          se = list(mod.1.a.se,mod.1.a.no, mod.1.a.se.bs, mod.1.b.se,mod.1.c.se,mod.1.d.se,mod.1.e.se),
          omit = c("country_id","year"),
          covariate.labels = c("Minority Cabinet","Kyoto Protocol",
                               "GDP Growth","Cabinet's Mean Environmental Protection",
                               "Total Greenhouse Gas Emissions/Capita","Green Party in Government",
                               "Green Party's Environmental Protection", "Post-Fukushima", 
                               "Share in Population worried about Environment",
                               "PM CIP net of Green Parties",  "Green Party's Polls",
                               "Green CIP gross CIP", 
                               "PM net CIP $\\times$ Green Party CIP","Intercept"),
          add.lines = list(c("Number of Countries","9","9","9", "9", "9", "9"),
                           c("Fixed-Effects", "Yes", "Yes","Yes", "Yes", "Yes", "Yes", "Yes"),
                           c("Trend", "Yes", "Yes","Yes", "Yes", "Yes",  "Yes","Yes"),
                           c("SEs", "Clustered", "Not Clustered", "Bootstrapped",  "Clustered", "Clustered", "Clustered", "Clustered")), dep.var.labels = "")


